Sediment delivery to sustain the Ganges-Brahmaputra delta under climate change and anthropogenic impacts

The principal nature-based solution for offsetting relative sea-level rise in the Ganges-Brahmaputra delta is the unabated delivery, dispersal, and deposition of the rivers’ ~1 billion-tonne annual sediment load. Recent hydrological transport modeling suggests that strengthening monsoon precipitation in the 21st century could increase this sediment delivery 34-60%; yet other studies demonstrate that sediment could decline 15-80% if planned dams and river diversions are fully implemented. We validate these modeled ranges by developing a comprehensive field-based sediment budget that quantifies the supply of Ganges-Brahmaputra river sediment under varying Holocene climate conditions. Our data reveal natural responses in sediment supply comparable to previously modeled results and suggest that increased sediment delivery may be capable of offsetting accelerated sea-level rise. This prospect for a naturally sustained Ganges-Brahmaputra delta presents possibilities beyond the dystopian future often posed for this system, but the implementation of currently proposed dams and diversions would preclude such opportunities.

. Grain size and bulk strontium data from Transect G sediment cores across the lower Ganges-Brahmaputra delta. Location shown in Figure 1. Upper panel plots the bulksediment strontium (Sr) concentration (this study), and the lower panel shows bulk sand-mud sediment distribution and calibrated radiocarbon ages (Grall et al., 2018) 1 . Also shown are Pleistocene surfaces that define the three lowstand river valleys. Note that Sr data readily distinguishes the provenance of Ganges, Brahmaputra, and Meghna sourced river sediments. Strontium values reveal the stratigraphic units of the three rivers to be largely discrete with limited mixing, enabling us to calculate the total mass of delta sediments contributed by each river.

"Time-equivalent Units"
Mapping of the "time-equivalent units" (TEQ) was performed using ArcGIS. Published maps of subsidence 1 (Fig. S2) and eustatic sea level 2 were combined to reconstruct the accommodation produced by relative sea-level rise and available for infilling with delivered sediment. Extensive prior research has shown that the sequence of deltaic facies in the G-B is dominated by aggradation (not translation), requiring that accommodation was readily infilled during Holocene sea-level rise and thus can be used to faithfully reconstruct the spatiotemporal distribution of sediment storage and serve as a proxy for fluvial sediment delivery 1,3-16 . We then calculated the areas and volumes of overlap between the Pleistocene-Holocene boundary and the reconstructed surface of relative sea-level change (i.e., TEQ) for each time unit (12-10 ka, 10-8 ka, 8-6 ka, and 6-0 ka) (Fig. S2). The Holocene-Pleistocene boundary depth was first mapped by hand using grain size, oxidation, and radiocarbon data from the BanglaPIRE dataset 1 and supplemented with additional core data from many other studies 3,9,10,[17][18][19][20][21][22][23][24][25][26][27] . Simple cokriging was used to create an interpolated prediction surface of the depth to the Pleistocene boundary ( Fig. S3) from which the total volume of the Holocene sediment package was calculated. We then subtracted each 3D TEQ surface from the Holocene-Pleistocene boundary surface to calculate the volume of sediment stored on the delta within each time period (Table  S1). We isolated the cores and samples within each TEQ to estimate the grain-size distribution of the sediment and converted sediment volume to mass using measured bulk densities (1.   Holocene sediment sequences are thickest in the river valley and basin settings (like Sylhet basin, NE corner of the delta), but the Pleistocene surface is only shallowly buried elsewhere on the delta, particularly near the delta apex. The red line shows the boundary between the fluvial delta (slope 10 -4 ) and the fluvial-tidal delta (10 -5 ) 12 . The World Ocean Base layer is from Esri, GEBCO, NOAA NGDC, HERE, Garmin, and other contributors. Table S1. Volume of sediment stored on the delta within each time period (TEQ). Bulk densities are relatively higher from 10-8 ka and 12-10 ka due to deep burial and a primarily sandy composition. Bulk densities for sediments deposited between 8-6 ka and 6-0 ka are relatively lower due to their shallower burial depths and to greater fraction of mass stored as lessconsolidated muds in the offshore subaqueous delta. 12-10 ka 909 1800 1636200 818 * Offshore volume distributed with 75% of offshore material contained in the 6-0 ka unit and 25% of offshore material contained in the 8-6 ka unit. Bulk density applied to offshore material is 1,100 Mt/km 3 .

Reconstructed potential sediment loads
Reconstructed suspended and bedload estimates of the Holocene sediment load using the observed ratio of mud-sized sediments to sand-sized sediments in shallow (<~30 m) sediment (last 6 ka -see methods) and from Akter et al. (2016) 11 to account for preservation biases are shown in Tables S2 and S3. Thus, total fluvial sediment load during the early Holocene may have been 30-70% higher than the minimum storage rates reconstructed here, thus potentially increasing sediment delivery estimates by ~600 Mt/yr for 12-10 ka and 10-8 ka. When sediment loads are reconstructed using the suspended:bedload data from Akter et al. (2016) 11 , the overall Holocene load would be ~5% higher but the last 6 ka would decrease ~25%, suggesting a preservation bias for suspended load in these younger sediments, which have not yet been buried below the depth of channel reworking. The bedload fraction drops to 45% of the stored load between 12-10 ka and 30% between 10-8 ka to ~25% between 8 ka to present. These results suggest that the Ganges-Brahmaputra delta (G-B) fluvial sediment load may be coarser than previously recognized with an average of 20-25% bedload over the Holocene. Moreover, these estimates of a coarser and larger sediment load align with sediment stored offshore (on the order of several Mt/yr), suggesting that sediment discharge rates to the G-B were likely higher than those recorded by storage alone. Table S2. Potential sediment loads reconstructed using the ratio of total sediment (excluding bedload) and bed material estimated from Akter et al. (2016) 11 . Applying the average ratio of suspended sediment to bedload to measured bedload stored in the Holocene delta increases the potential sediment delivery over the last 12 kyr, particularly in the early Holocene (12-8 ka). Table S3. Potential sediment loads reconstructed using the ratio of mud-sized sediments to sands observed in sediments stored on the delta between 6-0 ka. Using this method, sediment loads could potentially exceed stored amounts by 20-50+%.

Error analysis and uncertainties in estimates
There is no tractable process for a formal error analysis on the study's results given the many steps from collecting core samples to converting the raw data to a finished mass distribution. The propagation of errors would hold little value for assessing the accuracy of results. However, to control for uncertainty through the budget analysis, we used results from the physiographic regions with the highest data density as controls for areas with less data, such as West Bengal. These comparisons generally confirmed that the initial results were appropriate, and the few adjustments made were relatively minor and are noted in the main text. We compared our observations with those from the many studies cited in the paper, and we incorporated their data into the budget calculations wherever possible. We also evaluated uncertainty by varying key values such as the bulk density conversions and distribution of the Holocene boundary. These semi-quantitative assessments typically resulted in budget and mass deviations of 5-15%.
Other known errors include small uncertainties due to clipping and other data processing in ArcGIS, although these are minor since they occur at the boundaries of the basin where Holocene sediment storage is limited by the shallow Pleistocene surface.
Although the budget may lack some precision, the overall accuracy of results appears robust. For example, we have high confidence in the major boundaries of the Holocene delta, which are well defined by highly oxidized paleosols, lowstand gravel surfaces, radiocarbon ages, and sharp upland transitions. These robust constraints mean that our total sediment calculations and mass delivery rates are well supported. In all, the sources of error and uncertainty in the budget prove to be relatively small compared to the spatial and temporal scales at which the system is being analyzed here. Based on our semi-quantitative assessments of uncertainty, we suggest that the error for any given value presented in the paper is ±10% at the 1-sigma level and ±20% at the 2-sigma level.
The greatest uncertainty may lie with the mass of sediment that bypassed the delta and was deposited either on the deep-sea Bengal Fan or advected laterally along the shelf. We evaluate this potential error in the reconstructed sediment loads presented in Tables S2 and S3, and we discuss these issues thoroughly in the main text. The potential sediment bypassing is most likely to affect values for the 12-10 ka and 10-8 ka TEQ periods, when sea level was lower and the rivers were constrained to their lowstand valleys. After 8 ka, though, there is little evidence that significant bypassing occurred once the rise of sea level slowed. Bypassing would also primarily affect the total mud mass, since fine-grained sediments will have greater transport distance and are more likely to be reworked within the delta (i.e., preservation bias). The high mud:sand ratio for the 0-6 ka and 6-8 ka, however, suggests that most of the finegrained sediments are preserved in the delta for these periods. We also restate the key point